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INTRODUCTION 


The objective of this research is to develop efficient numerical techniques for 
the study of aeroelastic response of a propfan in an unsteady transonic flow. A 
three dimensional unsteady Euler solver, developed at Georgia Institute of 
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Application of an Efficient Hybrid Scheme for Aeroelastic Analysis of 

Advanced Propellers 

Abstract 

An efficient three dimensional hybrid scheme is applied for solving Euler equa- 
tions to analyze advanced propellers. The scheme treats the spanwise direction semi 
explicitly and the other two directions implicitly, without affecting the accuracy, as 
compared to a fully implicit scheme. This leads to a reduction in computer time 
and memory requirement. 

The calculated power coefficients for two advanced propellers, SR3 and SR<L, 
and various advance ratios showed good correlation with experiment. Spanwise dis- 
tribution of elemental power coefficient and steady pressure coefficient differences 
also showed good agreement with experiment. A study of the effect of structural 
febdbility on the performance of the advanced propellers showed that structural 
deformation due to centrifugal and aero loading should be included for better cor- 
relation. 


Introduction 

It has been known for some time now that the best propulsive efficiency is offered 
by propellers. However the efficiency drops off very rapidly as the cruise Mach 
number increases beyond 0.5, as high tip Mach numbers lead to high compressibility 
losses (due to wave drag). Currently an effort is underway to improve the propulsive 
efficiency of commercial and military aircraft. Newly designed high speed advanced 
propellers, also known as propfans, show a very high propulsive efficiency at cruise 

speeds upto Mach 0.8 [1]. 

The propfans are designed to delay the compressibility losses, thus extending 
the high efficiency of a propeller to relatively higher cruise Mach numbers. This 
is accomplished by sweeping the blade backwards and using thinner airfoils, on 


1 


improved on this method by using the curved lifting line concept to account for 
the sweep. In this approach the vortex wake is represented by a finite number of 
vortex filaments in place of the continuous sheet of vorticity as used in Goldstein s 
approach. The analysis has been further extended in reference [5] by placing the 
vortex filaments along the stream surfaces so that they conform to the shape of the 
axisymmetric nacelle. 

Hanson [6] and Williams [7] applied the Kernel function approach to a propfan 
blade. They numerically solve a linear integral equation for upwash angle due to the 
blade pressure distribution by discretizing the load representation. The friction drag 
is obtained from the two- dimensional airfoil tables as a function of lift coefficient for 
the appropriate section camber, thickness and a Mach number adjusted for sweep 
and three-dimensional effects. The induced drag is obtained by determining the 
kinetic energy-per-unit-length of the far wake. The methods mentioned so far are 
based on linearized analyses. However, as the advanced propeller operates at or 
near transonic tip Mach number, flow nonlinearties may become important. 

Jou [8] has applied the finite volume approach of Jameson [9] for the analysis 
of propfans using full potential equation. The formulation was not able to provide 
converged solutions for free stream Mach numbers greater than 0.6. It was concluded 
that strong rotational flow effects were present near the leading edge, which could 
not be modelled by the potential equation. In addition the potential flow equations 
at times, lead to non-unique solutions. 

Chausee [10] and Whitfield et al. [11] have applied the unsteady, three dimen- 
sional Euler equations to the propfan geometry. Matsuo et al. [12] have recently 
solved the full Navier - Stokes equations around a propfan. Some of these methods 
have been reviewed in reference [13], with regards to performance prediction. 

All the analyses mentioned so far, with the exception of Whitfield et al [11] 
have been for axisymmetric flows. For a propfan in flight configuration, the flow is 
not axisymmetric. Even for cruise conditions the nacelle is at an angle of attack 
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explicitly requires only two costly inversions of block tridiagonal matrix, as opposed 
to three inversions for a fully implicit scheme, per time step. It also reduces the 
memory requirement as only two time levels of information needs to be stored at 
any given time, one of which needs to be only two dimensional. The use of such 
a hybrid scheme leading to reduction in computer time and memory requirement, 
makes the scheme more efficient. 

The specific objectives of the present paper are 1) to apply an efficient hybrid 
scheme to analyze advanced propellers, 2) to calculate steady performance, 3) to 
include structural deformation, due to centrifugal and steady aero loading in the 
analysis, 4) to study the effects of structural flexibility on the performance of ad- 
vanced propellers. The governing equations and the numerical solution method are 
described first followed by results and discussion. The methods developed here are 
expected to be helpful for future aeroelastic research. 

Formulation 


Aerodynamic Model: 

The Euler equations, in conservation form, in Cartesian coordinate system can be 
written as: 


(q) f + (E), + (F) w + (G), = 0 


( 1 ) 


where q is the vector containing conserved flow properties. E,F and G are the 
nonlinear flux vectors which are functions of the vector q. The subscripts denote 
the partial derivative of the vector. In the above equation 
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( 4 ) 


rj = T}(x,y,z,t) 

C = £(*»»»*>*) 

T — t 

These coordinates are non orthogonal and completely general. The equation 
(1) can be rewritten as: 


i T + E € + F„ + G c - 0 


( 5 ) 


where 
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U,v, and W are the contravariant velocities, and J is the jacobian and £*, Tjx, 
<* etc. axe the metrics of transformation. 

Initial and Boundary Conditions 

A large number of problems can be described by the same set of governing 

Pfcwi-'EEtsf'iv-* DLAN’.\ MOi r 
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inflow boundary, all quantities are fixed to that of the free stream, as disturbances 
cannot travel upstream in a supersonic flow. At the subsonic outflow boundary, four 
characteristics should escape, thus the four quantities p, pu, pv, pw are extrapolated 
from inside while the pressure is fixed to that of the free stream. For supersonic 
outflow, all characteristics should escape, hence all quantities are extrapolated from 
inside the flow domain. 

The block interface boundary : 

It is neither efficient nor practical to solve all the blades at the same time, 
hence, one blade passage is handled at a time. This introduces additional bound- 
aries for computation. Across these boundaries all the variables must be continuous, 
except on solid boundaries and boundaries downstream of the blade. The bound- 
ary condition, for these boundaries, depends on the type of flow being solved. An 
axisymmetric flow would require periodicity on the fluid interface boundaries. Peri- 
odicity will require that the two boundaries have same fluid properties. As shown in 
figure (la), the fluid properties at the boundaries K=1 and K=KMAX are updated 
as the average of fluid properties at K=2 and K=KMAX-1 for a symmetric flow 
field. 

For an unsymmetric flow, the periodicity on these boundaries does not exist. 
Therefore, in order to obtain the solution for such a case, the whole propfan should 
be solved. This is done by advancing the solution of each block one time step, 
one block at a time. In this case again the boundaries are updated explicitly, after 
the interior points have been updated.This is done by averaging the flow variables 
from the nodes on each side of the boundary from the adjoining blocks. Referring 
to figure (lb), (the subscripts refer to the corresponding block) the quantities at 
boundary K=KMAX of block N (which is also the boundary K=1 for block N H- 1 ) 
would be the average of flow quantities at K=KMAX-1 of block N and K=2 of block 


.T'N: 
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marchine direction is reversed after every sweep, in order to remove any dependency 
on the marching direction. Equation (10) can then be rewritten as : 

q n+i _ q" _ At (E” +1 + F" ,T1+1 + G" +1 ) (H) 

Since the ry marching direction is changed every iteration, the F"’" +1 alternates 
between 


F n _ P" + l , 

ij + l.fc 

2Ary 


during the odd time steps, and 


2Ary 

during the even time steps. 

The above discretization leads to a set of algebraic equations for q. These 
equations are costly to solve since the flux vectors E and G are highly nonlinear. 
The nonlinearity is removed by linearising the fluxes about the previous time step 
value, resulting in the following linear equation : 


[J + At( 6 4 A" +6 c B n )]q n+1 = [/ + Ar(6^ n + 6<B n )]q n + R n,n+1 U 2 ) 

where 

R n,n+1 = _ At -f 6„F n ’ n+1 + 6 c G n ) (13) 

and the operator notation ^(Aq) = [^A]q and 6 c (Bq) = [6 C B ] q is used. 

This Euler equation formulation can be very easily extended to solve the Navier- 
Stokes equations by simply adding the viscous terms to the right hand side. This 
does not alter the numerical formulation. 
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name ‘Alternating Direction’. These inversions are performed at each spanwise 
station, marching along the spanwise direction. As mentioned earlier, the marching 
direction is reversed every iteration. Each element of the block tridiagonal matrix 
has 5x5 elements. 

Artificial Dissipation: 

The use of central difference, makes the scheme mildly unstable, and also introduces 
odd even decoupling. This is remedied by adding artificial dissipation. The imple- 
mentation of artificial dissipation, in the present work is based on the formulation 
of Jameson et al. [16]. This scheme has a second order implicit dissipation and a 
blend of second/ fourth difference explicit dissipation terms. A scaling factor for 
both implicit and explicit dissipation is employed to control the amount of dissipa- 
tion in the scheme. Adding the dissipation terms, equations (18) and (19) can be 
written as: 


[/ + Ar ( 6 ( A n + £/ D; { )] Aq* n+1 = R n ’" +1 - t E D E Ar (20) 

[/ + Ar(6 < B n + 6 / D / .)] Aq" +1 = Aq* n+1 (21) 

where D/ { and Dj ( are second order implicit dissipation terms and D E is the explicit 
dissipation term, given in reference [22]. ej and t E are user supplied constants, which 
depend on grid spacing. At the boundaries the fourth order differences are repalced 
by second order differences. 

AerQ< ? la§tk..MQdgl ; 

As mentioned earlier, the propfan has thin, swept, and twisted blades. Since the 
blades axe thin and flexible, deflections due to centrifugal and steady aero loads are 
large. Hence, the aeroelastic problem is inherently nonlinear, requiring geometric 
nonlinear theory of elasticity [17]. 

The blades have large sweep and twist, which couples blade bending and tor- 
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i.e., until the change in deflection from the ( i + l) th iteration is equal to that from 
the i th iteration. 


Results and Discussion 

The hybrid numerical scheme discussed in the previous section, was first applied 
to an isolated aircraft wing in reference [21] and to a helicopter rotor blade in 
reference [22]. Typical results showing blade loading, are reproduced in figures (3) 
and (4). As can be seen from both these figures, the hybrid scheme is able to predict 
flow phenomena of varying complexity with fairly good degree of accuracy. 

The propfan blade has a much more complex shape than the aircraft wing or 
the helicopter blade. The high twist, large sweep, low aspect ratio, close proximity 
of other blades, presence of nacelle and thinner blades near the tip, make the flow 
field around it very complex. In the following, the flow solutions obtained for 
two advanced propellers, namely SR3 and SR7L, are presented. The calculations 
have been performed on a ’hot shape’, obtained by including the deflections due to 
centrifugal loading in the undeflected blade shape (’cold shape’). 

A body fitted H-0 grid was used for these calculations. A typical grid used in 
the calculation is shown in figure (5). The domain of calculation was taken to be the 
region between two blades with upper surface of one blade and lower surface of the 
adjoining blade as the boundaries of the domain. This region is referred to as blade 
passage. In general, in order to model the influence of adjacent blades (cascade 
effect) the entire propfan with all the blades (passages) are solved. However, for an 
axisymmetric flow field, considered here, all blade passages can be assumed to be 
identical, and only one blade passage is solved enforcing the conditions of symmetry. 

SR3 Propfan 

The hybrid scheme, described earlier, was used to solve the flow field around 
an 8-bladed SR3 propfan. The SR3 propfan was designed to operate at a free 
stream Mach number of 0.8, advance ratio of of 3.06, at an altitude of 30,000 feet. 
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an overprediction in the tip region 


results in an under prediction on the inboard 


region. 

SR7L Propfan 

The SR7L propfan has been designed for an operating free stream Mach number 
of 0.8, rotational speed of 1700 rpm, at an altitude of 35,000 feet. In this section 
calculations for a two bladed SR7L propfan are presented. The aerodynamic calcu- 
lations are first performed on the ’hot shape’. The effect of blade flexibility is then 
included in the calculations. 

In figure (8) the elemental pressure coefficient difference is compared with ex- 
periment for a 2-bladed SR7L propfan. The blade was operating at a free stream 
Mach number of 0.775 and advance ratio of 3.088. The 75% span setting angle 
was adjusted to match the power coefficient by a rigid body rotation of the blade 
about the pitch change axis. The pressure coefficient difference A C p (C pt - C P J 
is plotted and compared against experimental data [27] at various span locations. 
The comparison is good, except near the leading edge on the outboard stations. 

The effect of blade flexibility on performance was studied next for the SR i L 
propfan blade. The effect of flexibility is included by the aeroelastic iteration pro- 
cess, described earlier and shown in figure (2). 

It is important that the blade finite element model accurately reflects its struc- 
tural characteristics, since the entire analysis process is centered around the stiffness 
matrix. The NASTRAN finite element model used in this study is based on the final 
blade design [24]. The SR7L blade has an aluminum spar, nickel sheath, and fiber 
glass shell with foam fill. The shell, adhesive, spar, and shell filler material were 
combined using the Composite Blade Structural Analysis (COBSTRAN) program 
to produce equivalent, monolithic shell elements [25]. The finite element model of 
the SR7L blade is shown in figure (9a). The model has 261 nodes and 449 triangular 
shell elements. Bar elements are used to model the shank. Multipoint constraint 
grid chords are used to define the shank/blade interface [26]. 
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the final blade shape. 

In figure (11) the thrust coefficient is plotted against the power coefficient for 
subsequent iterations. The setting angle used in the calculations has been obtained 
by rigid body rotation of the hot shape so as to match the power coefficient obtained 
by experiments. The experimental point is also plotted. The power coefficient ob- 
tained from the hot shape (marked 1), compares well with the experiment. However, 
the power coefficient changes considerably (marked 2 to 4), as the blade is allowed 
to deform under this load. It can be seen from figures (10) and (11), that the initial 
change in shape, lead to large change in power coefficient. For this particular case, 
under which the blades are loaded heavily almost 40% change in power coefficient is 
observed when the effect of aerodynamic loading is included in the blade shape. The 
subsequent changes are not as large. Hence, in order to obtain a better comparison 
with experimental power coefficient and load distribution on blade, the blade setting 
angle should be chosen such that the converged shape power coefficient is compared 
against the experimental data. This requires some trial and error m selecting the 
’cold’ or ’hot’ shape from which the aeroelastic iteration should be started. Arriving 
at the final blade shape might be critical for vibration and flutter calculations, as 
well. 

In figure(12) the relative change in twist angle over the span is plotted. This 
shows that the largest deflection occurs near the tip, with practically no deflection 
on the root sections. Also it should be noted that the variation in the blade twist is 
nonlinear and is largest near the tip. A rigid body rotation of the blade to account 
for the change in twist, would result in a linear variation along the span. This 
clearly shows that for better performance calculations, structural flexibility should 
be included in the analysis. 

Figure (13) shows the in- plane deflection of the blade planform and figure (14) 
shows the out of plane deflection of the blade at constant chord. Again, the largest 
deflection is towards the tip, with practically no deflection towards the root. Figure 
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4. The effect of aero loads was to compensate for the untwisting due to centrifugal 
loads 
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Figure-5a H-grid in Streamwise Plane 
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Figure-5b O-grid in Azimuthal Plane 
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Figure-7 
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Figure-9a Propfan Blade Finite Element Model Figure-9b SR7L Frequencies - vs. - RPM 
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Figure- 12 Change in Blade Setting Angle versus Blade Span 
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Figure-13 In Plane Deflection of Blade Planform 
Moo = 0.775, J = 3.088, 0 = 54.6 
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Figure- 14 Out of Plane Deflection at Mid Chord 
Moo = 0.775, J = 3.088, p = 54.6 





